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Abstract 

We explore the zero-temperature behavior of an assembly of bosons 
interacting through a zero-range, attractive potential. Because the two- 
body interaction admits a bound state, the many-body model is best 
described by a Hamiltonian that includes the coupling between atomic 
and molecular components. Due to the positive scattering length, the low- 
density collection is expected to remain stable against collapse despite the 
attraction between particles. Although a variational many-body analysis 
indicates a collapsing solution with only a molecular component to its 
condensate at low density, the expected atomic condensate solution can 
be obtained if the chemical potential is allowed to be complex valued. 
In addition to revealing two discrete eigenfrequencics associated with the 
molecular case, an expansion in small oscillations quantifies the imaginary 
part of the chemical potential as proportional to a coherent decay rate of 
the atomic condensate into a continuum of collective phonon excitations 
about the collapsing lower state. 
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1 Introduction 

This paper is intended to provide a complete quantum mechanical descrip- 
tion of the instability inherent in zero-temperature atom-molecule Bose sys- 
tems with attractive interparticle interactions, but with positive scattering 
length. At low-density it is well-known that the pressure in these systems 
is proportional to ap 2 , where a is the s-wave scattering length and p is the 
number density. This suggests that despite being innately attractive, the 
ensemble can remain stable against collapse if there is a nonnegative scat- 
tering length [1]. Nevertheless, we find a ground state that tends toward 
collapse, with the expected, positive-pressure case realized only at the cost 
of including an inherent quantum instability [2, 3, 4]. Elucidation of this 
result is provided by a series of derivations presented in the following sec- 
tions. 

Because many-body interactions are built up from the pairwise sum of 
those between individual particles, Section 2 begins with a review of the 
low-energy two-body scattering. After relating the scattering length to the 
interaction strength, it is shown that a separable interaction not only cap- 
tures the low-energy physics, but yields analytically tractable results. In 
addition, we include coupling to a molecular Feshbach state, thus enabling 
a tuning of the effective scattering length through a resonance. Solving the 
two-body problem provides a convenient context to model the many-body 
atom-molecule system, which is the subject of Section 3. As discussed in 
Section 4, a Gaussian variational procedure is used to find the lowest en- 
ergy states of the system. Upon deriving general variational expressions for 
these states, the static, uniform case is considered in Section 5 where it is 
found that a collapsing two-piece ground state persists even in the case of a 
positive scattering length. 

Notwithstanding, the expected solution can be obtained by allowing the 
chemical potential to assume complex values, where the imaginary part 
quantifies the decay of the condensate. Section 6 demonstrates this result 
for the special case of a uniform system whereas Section 7 proves the persis- 
tence of this feature even when the uniformity constraint is removed. In the 
final section, we expand in small oscillations, obtaining a pair of discrete ex- 
citation frequencies corresponding to the collapsing solution. Additionally, 
the expected energy per particle resides within a continuum of excitations 
above the collapsing ground state. Energy conservation indicates that the 
decay represents a transition into the continuum, thus revealing that the 
expected solution evolves into states that tend toward collapse as well. This 
analysis completes the result as it reveals the imaginary part of the chemi- 
cal potential to be associated with a decay rate into collective modes of the 
collapsing lower state. 
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2 Two-Body Scattering 

Because microscopic models of many-particle assemblies are built from all 
pairwise interactions between the constituents, a reasonable place to begin is 
by considering the interaction just between two particles. After identifying 
the s-wave scattering length as the relevant low-energy scattering parameter, 
this quantity is then related to the interaction strength. In doing so, we 
recognize relationships that any model potential must satisfy in order to 
accurately describe the interparticle interaction. From this analysis, we 
find that the separable potential is simple enough to be of calculational 
advantage, yet nonetheless retains sufficient flexibility to depict low-energy 
scattering. 

With the separable interaction, the scattering description is extended to 
include the possibility of molecular binding, a crucial element in the loss 
mechanisms at work in atomic condensates [5, 6]. Due to the presence of 
the hyperfine-induced molecular state, there exists a Feshbach resonance 
that occurs when the effective scattering length becomes unbounded upon 
tuning an externally applied magnetic field. In completing the analysis, 
an effective range equation is derived for this coupled system. From these 
two-body results, a many-body Hamiltonian may be postulated in which 
subsequent investigation reveals the properties of the collective. 

2.1 The Scattering Length 

In the center of mass, two-body scattering processes are equivalent to re- 
duced masses, m/2, impinging on a force center that is identical to the inter- 
particle interaction. Within the influence of this potential, V, the reduced 
mass obeys the full time independent Schrodinger equation, 



where the Hamiltonian is a sum of a free kinetic term plus the potential, 
H = Hq + V. Long before collision, the incoming particles lie far outside the 
potential's range, thus allowing the incident flux to be represented as plane 
waves obeying the free Schrodinger equation, 



It is convenient to relate the plane waves to the full scattering state, \tp), by 
introducing the T operator, 






(2) 




(3) 
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From this definition, the solution to the full Schrddinger equation (1) may 
be written 

|^ = |k) + G (2A; 2 + ^)T(2A; 2 )|k), (4) 

where the free-particle Green's function 1 has the coordinate space represen- 
tation 

r pik'-(x-x') 

(x|G (2fc 2 + ^)|x / ) = / 



k' 



2k 2 +ie- 2k' 2 

(5) 

2 gifc|x— x'| 
8-7T |x — x'| 

In this expression, the integral is evaluated with complex contour integra- 
tion in the standard way. 2 

Obtaining the asymptotic behavior for large |x| requires the usual expan- 
sion |x — x'| = x(l — x • x'/x 2 + • • • ), where the leading term is kept in 
the denominator, but the first two terms are retained in the exponential in 
(5). Substitution of this expanded Green's function into (4) leads to the 
asymptotic scattering wave function, 

Vk(x) ► e lk ' x - J-(£:x|T(2A; 2 )|k> — . (6) 

Comparing this with the usual form ^k( x ) e* k x — A(6, (p)e lkx /x, it is seen 
that, up to an overall constant, the matrix element is simply the scattering 
amplitude: 

^(6?,^) = ^x|T(2fc 2 )|k>, (7) 

where the polar angles are those between x and k. 

Since we will ultimately be concerned with zero-temperature bosonic sys- 
tems, it is useful to examine the low-energy form of the scattering amplitude. 
Figure 1 shows the radial wave function, which at "zero energy" has a linear 
form where the s-wave scattering length is the axial intercept: 

u ,k(x) -> x - a. (8) 

Comparing uq^/x with the asymptotic form of the wave function, the scat- 
tering length is identified as the zero energy T matrix element, 

8vra = (k|T(2fc 2 )|k / )| k=k , =0 . (9) 



lr To connect with the many-body analysis, all energies are on the scale of h 2 /2m, 
where m is the actual particle mass, not the reduced mass. Hence, the energy eigenvalue 
is denoted as 2fc 2 . 

2 We use the notation f = (27r)~ 3 f d 3 k. 

k 
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u ok (x) y 




Figure 1. At points outside the influence of the potential, 
V(x), this schematic depicts the s-wave (£ = 0) radial wave 
function (solid line) and the corresponding "zero energy" lin- 
ear asymptote (dashed line). Part (a) shows the free forms 
of both. For a repulsive potential, (b) indicates a negative 
phase shift [compared with (a)], as the wave function is "re- 
pelled" by the potential. The point where the zero energy 
line intersects the radial axis is the scattering length, a, which 
must be non-negative for repulsion. Note that as the range 
of V{x) shrinks to zero, so does the scattering length. For 
an attractive potential, (c) shows that the intercept can land 
on the other side of the axis, indicating a negative scattering 
length. As the well depth increases, (d) shows that the wave 
function can get "pulled in" far enough to give a positive 
value for a. Thus, for attractive potentials, the scattering 
length may have either sign. 

2.2 Separable Potential 

A two-body model interaction is sought that is not only consistent with the 
low-energy scattering physics summarized by Fig. 1, but that leads to ana- 
lytically tractable results as well. In particular, it is required that there exist 
the possibility for either positive or negative scattering lengths in the case of 
attractive interactions. To address these considerations, we use a separable 
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potential characterized by the center of mass matrix element depending on 
the strength, A, in addition to the form factor, /(k): 3 

<k|y|k') = A/(k)/(k'). (io) 

Multiplication of (4) by V, then applying definition (3) results in the Lippmann- 
Schwinger equation, T = V + VG$T, from which it is possible to relate the 
strength with the scattering length. Using the momentum space represen- 
tation of the Green's function, (k"|G (2£; 2 + ie)\k"') = 5(k" - k"')/(2k 2 + 
ie — 2k" 2 ), obtains the T matrix as 

/(k)A/(k') 
/(k" 



<k|T(2^)|k') = J r AJ }*> ■ (11) 



1 + A 
k 



2k" 2 - 2k 2 - ie 



Equation (9) then relates the scattering length to the strength by 

111 

+ 7> (I 2 ) 



8irab g A b ' 

in which b is identified as the range through 1/b = j f(k) 2 /2k 2 . Here, the 

k 

scattering length carries a subscript to distinguish the background value 
from the case where a Feshbach state is included. 

Equation (12) is consistent with the intuitive results depicted in Fig. 1, 
thus demonstrating that the separable form captures the low-energy scat- 
tering physics. According to (12), an attractive interaction (A < 0) can 
yield either positive or negative scattering lengths, which is qualitatively 
consistent with Fig. 1 (c) and (d). Alternatively, a repulsive interaction 
(A > 0) only gives rise to non-negative scattering lengths, which is as re- 
quired. Henceforth, we confine our attention to attractive interactions in 
which the range is allowed to vanish for simplicity. In this limit, we take 
/(k) — ► 1, with A approaching zero from below in such a way that a& g re- 
mains finite. 



2.3 Coupled Channels Scattering 

Up to this point, all that has been considered pertains to the case of a 
single channel outcome in which two particles enter in a particular state or 
channel, scatter, then reemerge in the same channel. Since there also exists 
the possibility for bound states to arise, an analysis is undertaken in which 
the energetically favored channel of the scattering state is coupled to the 

3 The zero momentum form factor is arbitrary, since its value does not alter any of the 
calculations. For simplicity, we therefore take /(0) = 1. 
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energetically unfavorable molecular state. This coupled channels analysis 
forms the basis of the description of the Feshbach resonance, a phenomenon 
in which the effective scattering length becomes unbounded. 

As shown in Fig. 2, the difference in electronic spin configuration gives 
rise to a distinct potential for each of the scattered and bound states [7]. 
Furthermore, the schematic illustrates the molecular binding energy, E, as 
the energy difference between the bound state and two free atoms, relative 
to the molecular (closed channel) potential. Analogously, the detuning, e, is 
defined as the same, but relative to the scattering (open channel) potential. 
Since both potentials arise due to a difference in spin states, the detuning, 
and hence the binding energy, can be adjusted using the Zeeman interaction 
in the presence of an applied magnetic field, B. Therefore, the detuning 
varies linearly with the field, 

e^e + AgB, (13) 

where Ag is the difference in the g factors between atoms and molecules. 
To include the coupling between the two states, it is most convenient 
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Figure 2. Schematic showing the potential of the closed 
channel (solid curve) with that of the open channel (dash- 
dot curve). The binding energy E and the detuning e are 
shown relative to the free level of each case. 

to define projection operators P and Q onto the respective scattering and 
bound Hilbert spaces, Jtfs an d [8]. For notational convenience, the 
appropriate projections of the full wave function, \ip), are given by 

p\il>) = |V P > 

(14) 

Q\if>) = \^ Q ). 
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Likewise, the various projections of the Hamiltonian are given by 
PHP = H PP QHQ = H QQ 

(15) 

PHQ = H PQ QHP = H QP . 

Using the standard relationships, P = P*, Q = and Q^Q + P^P = 1, the 
projection of the Schrodinger equation, H\ip) = E\ip), onto the scattering 
space gives 

H PP \ip p ) + H PQ \^ Q ) = E\ip p ). (16) 

Similarly, the bound space counterpart is 

H QQ \^) + H QP \ij p )=E\^). (17) 

Consider first the momentum space representation of the scattering space 
projection of the Schrodinger equation (16). In this case, the Hamiltonian, 
Hpp, contains both kinetic and potential pieces, which in momentum space 
are 

{k\H PP \iJj p ) = 2fc 2 V(k) + f (k|y|kXkO, (18) 

k' 

where f/>(k) = (k|^ p ). Consisting of a linear combination of the orthonormal 
basis functions, {\4> n ) G ^b}, the bound state projection is written 

\^)=J2cn\<f>n)- (19) 

n 

Thus, upon substitution of the separable potential, the momentum space 
representation of Eq. (16) is 

(2k 2 - E)i>Qt) + A/(k) f /(k')V(k') + Y^ c n^\HpQ\<Pn) = 0. (20) 

i 

Elimination of the coupling term requires an examination of the bound state 
projection. 

Since the {|</>n)} are eigenstates of Hqq, their eigenvalues are identified 
with detunings, e n , as in (13) and Fig. 2. Projecting the molecular equa- 
tion (17) onto one of the bound states, \4> n ), results in a solution for the 
coefficients, 

c n = --£±2 J /OW). (21) 
k' 

with coupling constants, a n , defined through 

iT Q p|k'} = J>a,/(k')|0,). (22) 
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Because it must vanish as the two-body strength is taken to zero, the cou- 
pling is represented by the product Xa n . Substitution of the coefficients into 
Eq. (20) gives 

{2k 2 - e) m + (\-\ 2 Y J j^e) /(k) / = °> ^ 23 ) 

thus illustrating the coupling to several bound states. Assuming the effect 
of the highest lying molecular state to be dominant over the ones below, a 
single bound state (denoted with the subscript 1) suffices for the description. 
Comparison with the usual single channel result (obtained when a\ = 0) 
reveals that the molecular state can simply be included by a shift in the 
strength, 

After setting the energy to zero, a combination of (24) and (12) arrives at 
an expression for the effective scattering length due to the presence of the 
molecular state, 

1 -± + J2>L,. (25) 



a ( e ) a bg e — Acq 
From the detuning (13), the scattering length can be put into the standard 
form, explicitly dependent upon the magnetic field, 

a(B) = a b g ( 1 - ^|j> (26) 

where the resonance width, AB, and the resonant field, Bq, are defined by 

AB = 8ira bg a 2 /Ag 

(27) 

B = (Xaf - e)/Ag - AB. 

In addition, the full T matrix is obtained by the substitution of (24) into 
(11), thus giving 



1+ A 



A 2 af \ f f(k" 



<k|T(2fc 2 )|k'> = > 2 2 I r J fn jk- 2 • (28) 



2k 2 J J 2k" 2 -2k 2 -ie 
k" 

Since a pole in the T matrix indicates a bound state, the binding energies 
are given by the values of k 2 where the denominator vanishes. After some 
algebraic manipulations, there emerges the following equation for the zeros: 

1 7o 2 / 1 /(k») 2 a\ = 

87ra bg 2 J k" 2 k" 2 + l 2 e + 2 7 2 - Xa 2 ' 1 J 

k" 
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with the binding energy defined as — k 2 = Jq. For negative background 
scattering lengths, the zero range (b — > 0) limit is taken in which A — » 0~ 
and /(k) — > 1. In this limit, (29) reduces to a cubic equation in the square 
root of the binding energy: 




By the integral in (29), only roots with 70 > are permissible. Figure 3 
shows the zero-range cases of the scattering length resonance along with the 
binding energy for the 85 Rb system. 




Figure 3. For the case of 85 Rb, (a) shows the full scattering 
length of Eq. (26), with the horizontal dashed line as the as- 
ymptotic background value and the solid vertical indicating 
the location of the resonance at 154.6 G. (b) With a vertical, 
dashed boundary line separating the region of positive from 
that of negative scattering length, a plot is shown of the 85 Rb 
binding energy corresponding to Eq. (30). Approaching the 
resonance, both the scattering length and molecular size be- 
come unbounded, indicating an increasingly weakly bound 
molecular state. Confirming this intuition, the binding en- 
ergy is seen to increase from zero on resonance. 



2.4 Effective Range 



Frequently, the partial wave analysis is employed where the determination 
of the phase shifts allows subsequent derivation of any other quantity of 
physical interest. In particular, the scattering length, a, and the effective 
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range, r e ff, are extracted from the expansion 

1 
2 



k cot J (fc) = -a' 1 + \r eff k 2 + ... , (31) 



where 5o(&) is the s-wave phase shift for particles of energy 2k 2 . Recall that 
the scattering amplitude is given by a sum over all partial waves, —A(9, if) = 
k- 1 J2e( 2 ^ + l)e^( fc ) smS e {k)P e {coB9). If only the s-waves contribute, this 
reduces to A = (87r) _1 T = — A; _1 e 4<5 °( fc ) sin So(k), thus suggesting that the 
effective range may be found by an expansion of the ratio of real to imaginary 
parts of T: 

cot«o(*) = ^^. (32) 

These parts are obtained from (28) after the integral is evaluated with 
Jdk'/(k 12 - k 2 ) = (2/fc)" 1 [ln |(fc6 - 47r 2 )/(fc6 + 4tt 2 )| + in]. Expanding the 
logarithms for small k, then taking the ratio in (32) obtains 



k cot 6o{k) = h 

a 

giving the effective range 



b -2 
167raf (e — Aa 2 ) 



2vr 3 



k 2 + ... , (33) 



r eff = -K- 3 b - 327ra 2 (e - Aa 2 ) . (34) 

Note that in the absence of the molecular state (qi = 0), the effective range is 
simply proportional to the range, b, which justifies the identification. In the 
zero range limit, the effective range becomes negative, r e ff — > — 32-7ra 2 /e 2 , 
due to the presence of the molecular state [9] . 



3 Many-Body Hamiltonian of the Coupled System 

Having reviewed the two-body scattering physics, we now discuss a corre- 
sponding Hamiltonian for the many-body case. In field operator language, 
this Hamiltonian is 

Aa (35) 

Here, the first term is the kinetic energy whereas the second is the sum over 
all pairwise contributions from the two-body interaction. Given by the mol- 
ecule number times the detuning, the next term accounts for the energy of 
the molecular condensate component. With a proportionality constant Aai 
and a molecular form factor F, the final term is the interaction energy due 
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to coupling of atoms to molecules and vice versa [10]. 

Before embarking on the analysis, it is worthwhile to attempt a simpli- 
fication of the above model. Suppose there were no two-body potential 
present. Having only a single Feshbach state with a coupling parameter g 
and detuning e, Eq. (23) becomes 

(2k 2 - E) m - ^/(k) I /(k')^(k') = 0. (36) 

k' 

For a detuning far enough away from the binding energy, |e| 3> \E\, the 
effective interaction strength is ~ —g 2 /^ which can be chosen to be equal 
to the original separable potential strength, A. Accordingly, an appropriate 
detuning, e = 2/ A 2 , and a coupling, g = ^—2/X, are defined such that 
the requirement |e| 3> \E\ is self-consistently fulfilled in the limit of zero- 
range attraction, A — > 0~. Replacing the two-body potential with a term 
that couples the atoms to an auxiliary molecular field, X, gives an effective 
Hamiltonian [11], 



H = £ *t?afl*p + e E + e E 

a/3 a a 

+ 7f E ^Afh^fl + 4= E XlK^^p + H. c. 

a/?7 a/37 



(37) 



Although we have replaced the two-body interaction by a single Feshbach 
state, it must be emphasized that the two Hamiltonians in (35) and (37) are 
in general different for nonzero range. 4 



4 Gaussian Variational Principle 

To obtain an equation of state for the coupled system, we employ a varia- 
tional procedure in which a Gaussian trial functional is used in calculating 
the expectation value of the Hamiltonian [12]. An extremization then de- 
termines the various solutions arising from the model. As a first step in this 
procedure, we decompose the field operators in terms of their corresponding 
coordinates, ip, <f>, x, and their conjugate momenta, tt, a), z>: 



*a 



-j= (i> a + m^j ; = -j= - iTTo^j (38a) 
-j= (j> a + iu a ^J ; & a = -j= (j> a - iu a ^J (38b) 



4 It turns out that that the two models are equivalent in the zero-range limit. 
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(38c) 



Where appropriate, only the atomic field is used in illustrating the required 
operator relationships since the molecular counterparts follow by a straight- 
forward comparison. As usual, the field operators' commutation relations 
are 



\Lrt \T>t 







(39a) 



which imply 



(39b) 



In the functional Schrodinger picture, the state vector, depends on 

the field V', 

|tf,t)->tt[^,t]. (40) 

Analogous to the coordinate space representation of single-particle opera- 
tors, the action of tp and 7r on the many-body state is 

r£ a |\M)-^tf[^,t] (41a) 



(41b) 



Explicitly, the coherent many-body state is an ansatz consisting of a product 
of Gaussian functionals, with one for each of the three fields: 



VW,<t>', X ',t} =iV^exp| [^^(t)-iS a/3 (t) 



a/3 



<%(t) 



+ »5> a (t)5<(t) 



(42) 



exp 
iV x exp 



E 
E 



5x' a (tY - iv a (t)5 X ' a (t) 
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where Nf, Ns and N x are the normalization constants and the fluctuations 
are given by the fields minus their mean values: 



Si// a (t) = i>' a 




8<j)' a {t) = 4>' a 


-Mt) 




~Xa{t) 



(43a) 
(43b) 
(43c) 

Unlike the molecular functionals, the atomic Gaussian has extra degrees 
of freedom as it is parametrized by a symmetric width composed of a real 
part, G~p/4 = Gp a /4, along with its canonical conjugate, — Yjp a . With 
their conjugates set to zero, the corresponding molecular widths are taken 
to be 1/2. This simplification is justified due to the absence of any and 
XX terms in the Hamiltonian. Under a variational analysis, general molec- 
ular widths simply reduce to the aforementioned values as anticipated in the 
ansatz (42). Because each term contains a balance of creation-destruction 
operator pairs, the Hamiltonian remains invariant under any phase trans- 
formation of the trial functional: 



I*,*) 



-i7V0(t)|^ a 



(44) 



Put another way, the invariance with respect to the phase angle, 6(t), must 
be connected with a fixed total particle number, N . Due to the continuous 
symmetry in 8, there must exist a mode of zero energy, otherwise known as 
the Goldstone mode. Upon conducting a small oscillation analysis in the 
random phase approximation (RPA), the presence of this zero frequency is 
explicitly demonstrated. 

In this formalism, the mean of any operator, O, is calculated from the 
functional integral 



(y,t\6\V,t) = J Vi;'^*[i(j',t}d^[^,t}. 



(45) 



For calculating the mean, it is useful to construct a basis of creation and de- 
struction operators for the atomic Gaussian in (42). As linear combinations 
of (p and 7r, these operators are found to be 



(46a) 
, (46b) 
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omitting the explicit time dependence for convenience. A direct application 
on the many-body state reveals that 



i 



V2 



c a (t)\$>,t) =0, 



(47a) 
(47b) 



thus verifying the construction. Using the commutators (39b) along with 
the symmetry of G and £ we have 



4(*).c/3(<)] =-\G aP {t). 



(48) 



Transforming to the {c^, c} basis, it is possible to easily calculate all quan- 
tities of interest. Inversion of (46) leads to the desired expression for the 
field operators 



**(*) = 2 £ 



C/3(t) 

(49) 



Employing Eqs. (47)-(49), all required mean values are calculated in the 
following: 

(*,t|#«(t)|*,t) = -^[^ a (t) + i7r a (t)] = * a (t) (50a) 

(tt, t\$ a (t)$p(t)\y, t) = -D af3 (t) + * a (t)*i3(t) (50b) 
t|*t (t)^(t)|*, t) = R aP {t) + (t)*^*) (50c) 

+ R^{t)R a& {t) + R ps {t)R ai {t) 

+ ^(^(iji^t) + ^(i)* 7 (t)ii a5 (t) (50d) 

+ **(t)^(t)^ 7 (t) + (tjtf^i^t) 

- ** (t)^(t)D 7 i(t) + #* (t)tt£(t)¥ 7 (t)tf,ft) 
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(^f,t\i—\^,t) = Yl [^a(t)^a(t)+U a (t)^ a (t)+V a (t)Xa(t) 



+ ^2 ^a0(t)Gi3a(t) +AfO(t) + total time derivatives. 

a/3 



(50e) 



In the expressions above, we have introduced R a p(t) and D a p{£) as the 
respective fluctuations of (^^) and (iff iff) about their mean field values of 
|\]/| 2 /2 and \E' 2 /2. Explicitly, these fluctuations are parametrized by G~^(t) 
and its conjugate, Y, a p(t): 



R a p(t) 
D a p(t) 



\G-%(t) + G a(3 {t)-8 a p +2^S Q7 (t)G 7(5 (t)S^(t), (51a) 
W a l(t)-G aP (t)\ + 2^£ (n (i)G 75 (t)£ 5/3 (i) 



7<5 



(51b) 



,J2 [^(t)G lP (t) + G a7 (t)S 7/3 (t)] 



Note that Eq. (50e) takes account of the wave functional phase, as given in 
the transformation (44). 

Lastly, Hamilton's equations of motion follow from the effective action 
which is defined as 



(52) 



S = J L(t)dt = J dt{V,t\idt-H\V,t). 
Substitution of the mean values (50) into the action obtains 

dt I [^a(t)^ a (t) + U a (t)4> a (t) + V a (t)Xa(t) 



a/3 J 



(53) 



where 



n = {y,t\H\v,t} 



(54) 



N = (^f,t\N\^f,t); 



(55) 
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a 

Since TV is constant, 6 must be time independent, thus giving 9 = =>■ # = 
const. = /i, which we identify as the chemical potential at zero temperature. 
By stationarizing the action, the equations of motion emerge as 

G a p[t) = 56a) 
£<*/?(*) = xr , 77s (56b) 

^ = fa^j (56c) 

faW = 7T77V— • ( 56d 

Not shown are the two molecular fields since their derivatives are straight- 
forward analogs of (56c) and (56d). From these it is evident that (M,9), 
(ir,ip), (w,4>), (f, X) an d (S, G) are canonical pairs. 

Employing Eqs. (50) and (55), the grand canonical Hamiltonian corre- 
sponding to (37) has the momentum space expectation value of 

K = H- i iM= j (k 2 - /i) [i?(k,k,i) + #*(k,f)*(k,t)] 
k 

+ (e-2/i) J $*(k,t)$(k,t) + (£-2 M ) y X*(k,t)X(k,t) 



k 



k,k',k" 

x { ~*(k", t) [-D(k, k', i) + (k, i)*(k', t)] 
+ E(k", t) [-D*(k, k', t) + #*(k, i)**(k', i)] } . (57) 

In the kinetic energy term, we have used the momentum space form T(k, k') = 
k 2 5(k — k'), whereas the coupling term's form factor, F, has been chosen 
to be proportional to / in the separable potential introduced in (10). Sim- 
plifying the notation, two molecular fields have been combined in the term 
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3(k", t) = Aai<3?(k", t) + gX(k", t). b Before trying to analyze the full time- 
dependent problem, it is beneficial to first examine the special case of the 
static, uniform medium. Moreover, the uniform results may be applied to 
nonuniform trapping geometries by an application of the local density ap- 
proximation. 



5 Static, Uniform Solution 

The static, uniform case has zero momenta (S = 0, ir = lo = v = 0) with 
constant mean fields in all of space, thus imparting a continuous translational 
symmetry to the system. For a hard sphere bose gas, the uniform energy 
per particle (on the scale of ti 2 /2m) is known to be 



e = ^=4™ P+ ^a 5 /y /2 + ... 
p 15 



(58) 



where u is the energy density and a is the s-wave scattering length [13]. It 
has also been shown that the lowest order term Airap is independent of the 
form of the two-body interaction [14]. 

We compare this result with a variational many-body analysis of the 
Hamiltonian in (57) which is further simplified by noting that for a uni- 
form system, the Gaussian trial functional's width is diagonal, 



G(k,k') = G(k)<5(k - k') 



(59) 



According to Eqs. (51), this yields diagonal forms for the fluctuations as 
well 



i?(k,k') 



iG(k)- 1 + G(k)-l 



<5(k - k') 



= i?(k)5(k - k') 



(60a) 



D(k,k') 



i G(k)" 1 - G(k) 



<5(k - k') 



D(k)5(k - k') 



(60b) 



$ and X are defined in analogy to (50a). 
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Likewise, in momentum space, the mean fields are simply constants multi- 
plied by delta functions: 

ip(k) = #(k) (61a) 
0(k) = <M(k) (61b) 

X(k)= x *(k). (61c) 

To obtain finite quantities, the Hamiltonian is divided by the volume of 
space. Applied to a uniform system, the First Law of thermodynamics, 
dE = —PdV + fidM, reveals the resulting quantity as the negative of the 
pressure, 6 

-P = i (« - «A0 = /(* 2 " C)«W + { J /(k)D(k) + (£ - (i) ^</> 2 

(62) 

1 Vi+aV+fi Vlv 



In addition to taking the mean fields as real quantities, we have defined 

£ > = \ ai (j) + g X (63a) 

j] = \ax<j> - gx- (63b) 

Along with G(k), the variational parameters include the mean fields tp, £ 
and rj. Extremizing on each results in the following set: 



- = => G(k) = - J ^7^7(k) ( 64a ) 



5P 



(f _ M ) ^ = o (64b) 

^ = o => ij^e + / /(k)D(k) + V = 0, (64c) 
or? r + a J 2 



where <r = (e - 2^)/(2Aai ) 2 and r = (e - 2p)/{2g) 2 . We use a step function 
form factor, /(k) = #(|k| — 47r 2 /6), so that both i?(k) and -D(k) vanish for 
> 4vr 2 /6 where G = 1/2. Thus, the upper limit in all radial integrals is 
cut off at 47r 2 /6 while taking /(k) = 1 everywhere. As a further reduction, 



6 The sign in front of -D(k) is unimportant since it can be absorbed into the as yet 
unknown parameter £. 
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the zero-range limit (b — > and A — ► ) is taken, whence we require the 
following quantities to the appropriate order in A: 7 



4tct 
r + a 

g + ?7 
2Xai 



1 



X 



2g 



A e - 2/i 
e-2/i 



(65a) 
(65b) 

(65c) 



After substitution of (64a) into Eqs. (60), the required integrals are writ- 
ten out explicitly, each being expressible in terms of elliptic integrals of the 
first and second kind, denoted as F and E, respectively [15]: 



k 2 R(k) 



1 

4^2 



47T 2 /fe 



fc 2 + 7 2 



^(£2 + 7 2)2 _ £2 



dk 



b^o 20vr 2 



[- (3£ 2 + 7 4 ) E + 7 2 ( 7 2 



(66a) 



F + 



2/) 



4n 2 /b 



R(k) 



Air 2 



k 2 + 1 2 



>4S S b 2 ^-(7 2 -lci)^] 



6^0 12vr 2 



(66b) 



4lT 2 /b 



b 



(66c) 



The resulting expressions have been expanded in the b — > limit, omitting 
terms of order 6 and higher. Identifying = — j 2 , the elliptic integral 
arguments are implicitly understood to be \/2|£|/(|£| + 7 2 ). 



7 Recall that g = sj-2j\ and e = 2/ A 2 . 
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Retaining only the nonvanishing terms in the A — > limit, the pressure 
and the density become 



i ate e 



-p = /(^-^(k) + e/^(k) + ( e-M)^ 2 



k 

e 

167ra(e — 2/i) 
J R(k) + \i? + 



(67) 



aft 



(68) 



As in (67), we will frequently find it convenient to use the full scattering 
length (25), where e is shifted by — 2/i. Both pressure and density depend 
on three parameters, tp, 7 2 = —p and £, which are constrained by the 
two variational Equations (64b) and (64c). Through the relation (68), the 
density is chosen as the free parameter. 



5.1 V = Solution 



First, note that Eq. (64b) admits two solutions, with one for ip = and 
the other for £ = p, = — 7 2 . We begin with the former, describing a system 
composed entirely of a molecular condensate and correlated atom pairs, but 
with no atomic condensate component. Substitution of (25), (65a) and (66c) 
into (64c) results in an expression that relates £ to 7 2 : 

4n 2 /b 



8vra(e + 2 7 2 ) 4vr 2 J 
o 

where, constrained by the form of the integral, £ ranges from — 7 2 to zero. 
Differentiation of the density expansion (68) indicates that dp/dC 2 > 0, 
thus implying an increasing density with increasing £ 2 . When £ = — 7 2 , 
a critical density is reached since (69) no longer admits a real solution for 
£ 2 > 7 4 . According to (68), dilute gases (pa 3 <C 1) in the low-density regime 
correspond to small £. In the limit £ — > 0, (69) reduces to 



dk = 0, 



(69) 



7oa(e + 2 7o 2 ) = 1. 



(70) 
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When compared with (30), this is identified as the equation for the binding 
energy 27^. 

In addition to (69), the equation of state for tp = is specified by the en- 
ergy and number densities which, when cast in terms of the elliptic integrals 
become 



-P + MP = m 2 ~1 A )E + 7 2 (7 2 - ^l) F] 

(71) 



16vra(e + 2 7 2 ) (e + 2 7 2 ) 2 

" = ^P^-^-l« f ] + (7Tl^- < 72 > 

From the limit £ — > 0, the energy per particle, e = u/p, reduces to minus 
half the molecular binding energy, —79- This result makes physical sense 
because at zero density, only the molecular state remains. Although the 
ip = solution terminates when £ = — 7 2 , the variational Equation (64b) also 
permits an alternative solution given by £ = fi = — 7 2 , which corresponds to 
a nonzero atomic field in general. 



5.2 V + Solution 

The solution for £ = — 7 2 permits an elementary evaluation of the integrals 
in (66): 

V^Tzt-l&t (73a) 



£>(k) > ^ - ^| 7 3 - (73c) 

6-»0 6 47T 2 



Unlike the tjj = case, (64c) is now solved for the atomic field: 

V — ^ + ^ 7 3 (74) 

2^ 6-0 8™(e + 2 7 2 ) ^ 4vr 27 { ] 
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With these expressions, both the energy and number densities are parametrized 



by 7 2 : 



u 



1 



16vra 69 2 ( e + 2 7 2 )' 



' = - s^W 2 + ^ 73 + (7TW 74 - (76) 

For the example of 85 Rb, Figure 4 shows the energy per particle curve 
formed by the merging of the ip = and i/j ^ solutions at the critical point. 
The fact that this point marks a quantum phase transition [16, 17, 18] 
follows from a discontinuity in the density between the two cases. That is, 
one piece has an atomic condensate density component and the other does 
not. Due to the persistently negative slope of the energy per particle, the 
pressure, P = p 2 de/dp, is always negative, by which we denote this two- 
piece solution as the "collapsing ground state" of the model Hamiltonian 
(37). Of particular importance is that this collapsing behavior persists even 
when the full scattering length is tuned to positive values, a regime thought 
to be stable against collapse. Indeed, this result challenges the intuitive low- 
density scattering length behavior, e ~ 4irap, indicated in expansion (58). 



6 A Complex Chemical Potential 

Due to an innately negative scattering length, 85 Rb has two-body interac- 
tions that are attractive, resulting in an ensemble that tends toward collapse 
as quantum degeneracy is approached. As such, this particular isotope can- 
not form a condensate containing more atoms than some critical number. 
However, by exploiting a Feshbach resonance, the appearance of the molecu- 
lar state allows the effective scattering length to be tuned to positive values, 
as indicated in Fig. 3(a). By tuning to positive values, a stable conden- 
sate has been experimentally observed in which the low-density energy per 
particle is found to have the usual behavior, e ~ 4nap, where a is the full 
scattering length [1]. 

Despite the positive scattering length, our previous solutions predict that 
a uniform assembly should still tend toward collapse, confronting us with 
the dilemma that the separable potential somehow does not capture the cor- 
rect many-body physics even though it does reproduce the correct two-body 
scattering (see Sec. 2.2). To resolve this issue, it is necessary to reexamine 
the solutions obtained in the collapsing case. Since the experiments have 
observed the formation of an atomic condensate, we focus on the solution for 
ip 0, exploring the consequences of extending this solution to zero density. 



26 QUANTUM INSTABILITY OF ATTRACTIVE BOSE SYSTEMS 



-1.5 




c -i./ 
^ BE/2 

I -1.9 

£ -2 




lu -2.4 



-2.5 



0.5 1 P c 1.5 2 2.5 3 
Density (cm -3 ) x10 16 



Figure 4. Figure showing the two-piece collapsing ground 
state for the case of 85 Rb at a magnetic field of 162.3 G, 
corresponding to a critical density, p c , of 1.22 x 10 16 cm -3 
with a full scattering length of 193 Bohr radii [26]. In the 
limit of zero density, the energy per particle simply reduces 
to half the molecular binding energy, BE/2 ~ 1.85neV. By 
tuning the magnetic field closer to resonance, the relatively 
high critical density may be brought within the regime of 
current experiments. 

To do so, the density Equation (76) is inverted, thereby obtaining fj, as an 
expansion in ^fp: 



When this expansion is substituted into (75), the usual low-density energy 
dependence, e = u/p ~ 4irap, is achieved. Although this does give the 
appropriate form, the chemical potential becomes complex at higher order. 
At first it may seem that this is simply an unphysical solution, but upon 
further consideration it is recognized as a signature of an instability that is 
inherent in the original system. 

With the Hamiltonian given by (37), consider the Heisenberg equation 
of motion for the atomic field operator 'ft: ihdt^l = [$\,H]. Employing 



the commutators [* a , = S a p, [^J,^^] = Stf $1 and [*J,^ 7 ^] = 
—Sfri&p — 5^/3 the right-hand side can be readily evaluated using the 




64a(e) 2 
3tt 2 



(77) 



+ 8yra(e)^- 64vr 2 a(e) 2 p 2 



GEORGE E. CRAGG AND ARTHUR K. KERMAN 



27 



expectation value of the mean atomic field (50a). This leads to 

d^* r- 
ih — = -y/2 (Aai$* + gX*) 

= -V2E*y. (78) 

We seek to add a phase to the mean fields such that (78) is consistent 
with the variational Equations (64) . If we let 

* = ^Ve" iM * /a , (79a) 

E = -^£e- 2i ^ h , (79b) 

then (78) gives back (64b). Thus, it is reasonable to interpret the chemical 
potential as the phase of the mean fields [19], where the imaginary part of 
H leads to a decay rate given by 

r = A^ 5 i2a(6) 5 /y /2 + ---- (80) 
2m 3 

Displaying unique dependencies on both the scattering length [~ a(e) 5 / 2 ] and 
the density (~ p 3 ^ 2 ), this coherent rate can be tested by further experiments. 
Under the conditions of the 85 Rb experiment, there were 10 atoms within 
a cloud of radius 25 fim when the scattering length was tuned to 193ao 
(ao = Bohr radius). These parameters yield a decay time, r ~ 1/T, of 14.3 
seconds, 8 which is in qualitative agreement with the observed 10 second 
lifetime [1]. 

Because all quantities depend on p, which is in general complex, it follows 
that the various other thermodynamic functions assume a complex character 
as well. This is not unphysical since the imaginary parts should simply be 
regarded as signatures of the coherent decay with the real parts assuming 
their usual physical interpretations. Although complex frequencies are well- 
known to signify the damping of collective modes [20], we emphasize that 
the distinctive effect here predicts a decay of the condensate itself, in the 
absence of any excitations. Nonetheless, a physical analogy may be found in 
the context of quantum electrodynamics, where the instability of the electric 
field emerges through a complex action [21]. 



For an attractive interaction, the longest decay time is attained in the limit A — > 
since the rate increases with decreasing A [4]. 
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Figure 5. For the case of Rb, we plot the real part of 
the energy per particle, e = u/p [see Eqs. (75) and (77)], 
corresponding to the decaying excited state (1). Below this 
is the two-piece collapsing ground state (2) as shown in Fig. 4. 
All results are for an applied magnetic field of 162. 3G. 

7 Decay in the Nonuniform Case 

The complex chemical potential persists even in the static, but nonuniform 
case. Thus far, only the uniform solutions have been discussed since it has 
been assumed that any nonuniformity can be accommodated using a local 
density approximation. Therefore, one could conjecture that a chemical po- 
tential assuming a complex value is simply an artifact arising because the 
uniform solution is too restrictive to capture the physics observed in exper- 
iment. 

To address this concern, the origin of the complex p must be identified. 
First, it is recognized that the expected solution has an equation of state 
associated with a positive chemical potential. At low density, the energy per 
particle and its derivative satisfy e ~ Anap > and de/dp > 0, respectively. 
By the relationship 9 p = du/dp and the definition of e as u/p, it follows 
that n = p{de/dp) + e. The positivity of e, p and de/dp thus implies the 
same for p. From (64a), it is seen that if p is positive, then G becomes 
imaginary for sufficiently small values of k. By construction, G" 1 is the real 
part of the width, raising an apparent contradiction unless the chemical po- 
tential is allowed to be complex instead. Therefore, an imaginary G is also 
a harbinger of the same instability described in the previous section. Most 
importantly, the ranging of G into imaginary values provides a convenient 



This can be seen by extremizing the pressure, P = /ip—u, with respect to p, 8P/Sp — 0. 
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test for this instability when discussing the nonuniform case where a full 
analytic expression for \i is not possible. 

For static nonuniformities, the grand canonical Hamiltonian can be writ- 
ten 



'a 



a/3 ^ ' ^ 'a 

^ ' a a/3"/ ^ ' 

(81) 

where h a p = T a p — fj, 8 a p. The fluctuation terms are given by the static form 
of (51): 

Rap = \ (jG-i + G Q/3 - sj\ (82a) 
D ^ = l(]G-l-G aP ). (82b) 



Differentiating with respect to G pa yields 
5R a(3 _ 1 / 1 



5G pcr 2 



iG^G^ + SapSfJ) (83a) 



8D n a 1 / 1 



a/3 _ 1 / J-„_i„_i 



5G pa 2 



J^G^-^A (83b) 



Using these derivatives, (81) is extremized on the width 

= ~ I E G p« ( h °0 + G- p l + IJ2 6 p° ( h "P - P "P) V = 0, (84) 

t" 7 afi a(3 

with F related to the form factor F by 10 

F a p = ^2Re{Cr}F T , a p. (85) 

After defining new matrices Y and Z as 

*ar/3 = K(3 + F a /3 (86a) 

= h a /3 — F a /3, (86b) 
the variational Equation (84) can be put into the compact matrix form 

--G^YG' 1 + Z = 0. (87) 
4 



3 We have used the symmetry F T - a p — F T - t p a 
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To obtain a symmetric solution for G 1 , we first multiply by VY on the 
right and on the left, 

--Vyg- 1 Vy Vyg^ 1 Vy + Vy z Vy = o. 

4 

Moving the first term to the right-hand side, then taking the square root, 
G^ 1 is seen to be 11 



G^ 1 = 2 —= y VY zVy —=. (89) 

Vy v Vy 

Likewise, the symmetric form of G is most easily obtained by multiplying 
(87) on the left and right by G, then following the same steps to end up 
with 

^i-Ly^f-L (90) 

In momentum space, the operator h is given by /i(k, k') = (k 2 — fi)5(k— k'). 
For k 2 = fj,, h = 0. At such values of k, Eqs. (86) show that Y = —Z = F, 
thus giving 



G 



fe 2 = 



F 



-F 



(91) 



2' 

indicating that G is in general complex. By our earlier observation, it follows 
that the coherent decay is present even in the nonuniform case. 



8 Small Oscillations: Generalized Random Phase Approx- 
imation (RPA) 

Evidence for the production of a coherent atom-molecule condensate has 
been demonstrated by the dynamical response of the 85 Rb system near a 
Feshbach resonance. In particular, oscillations in the atomic density were 
observed in trapped samples by Donley et al. [22], indicating the coexis- 
tence of a molecular counterpart. Although theoretical investigations have 
been conducted by Holland et al. [23] and by Kokkelmans and Holland 
[24], it is perplexing that only the latter seems to predict a much longer 
damping rate in agreement with the experiment. In any event, all of these 
investigations dealt with the dynamics of the expected solution which ex- 
pands against the trap due to the positive scattering length. Instead, our 

11 This series of manipulations obtains the required symmetric form for G" 1 . In general, 
however, this is not a unique solution since multiplying (87) on the right by Y, then taking 
the square root results in G _1 = 2 \/~ZY Y^ 1 , which is different than the solution in (89). 
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aim is to carry out a small oscillation expansion of the collapsing ground 
state since its spectrum completes the physical picture by providing an in- 
terpretation of the coherent decay uncovered in the previous two sections. 

To find the collapsing state's low-lying excitations, all variational quanti- 
ties are expanded about their stationary values, resulting in an oscillator-like 
Hamiltonian expressed in terms of effective mass (21) and spring (OS) matri- 
ces. Solving the full problem is accomplished in a series of stages, with each 
including an added generality over the previous. As the first and simplest 
step, we consider only the diagonal or noninteracting part of the Hamilton- 
ian from which the corresponding (5G, <5E) oscillations represent the energy 
of two free quasi-bosons. Including both the diagonal and off-diagonal el- 
ements in 21, a zero is obtained as a discrete point in its spectrum, thus 
verifying the presence of the Goldstone mode in general. Nevertheless, a 
complete accounting of the excitations must be obtained from the product 
21- 23. Correspondingly, the eigenvalue problem has the familiar Lippmann- 
Schwinger form, indicating the associated eigenfrequencies to be identified 
with the energy of two interacting quasi-bosons. In addition to discrete 
eigenfrequencies, there exists a phonon continuum of the collapsing solu- 
tion. Since the expected energy per particle lies within the continuum, en- 
ergy conservation gives a natural interpretation of the decay as a transition 
into the phonon excitations of the collapsing state. 



8.1 General Expansion 

All stationary quantities are expanded about their uniform equilibrium val- 
ues: 



G(k, k',i) 


= G(k)<5(k - k') + 5G(k, k',t) 


(92a) 


S(k,k',t) 


= <5E(k,k',t) 


(92b) 


V>(k,t) 


= #(k) + (ty(k,t) 


(92c) 


vr(k,t) 


= 5-K(k,t) 


(92d) 


(j)(k,t) 


= <M(k) + 5<p(k,t) 


(92e) 


u(k,t) 


= 5uj(k,t) 


(92f) 


x(M 


= X S(k)+5 X (k,t) 


(92g) 


z/(k,t) 


= 5u(k,t) 


(92h) 
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£(k,i) = tf(k) + 5£(k,i). (92i) 
It will prove convenient to define new momenta 

P = k - k' (93a) 

q=l(k + k'), (93b) 

having the interpretation as the respective total and relative momenta of a 
pair of quasiparticles. 

To illustrate the relationship between the variations and their inverses, 
consider the expression for G~ 1 (k, k',t), which, to second-order in 5G, we 
write as 

G- l (k,k',t) = G(k)~ 1 (J(k-k / ) + 5G~\k,k' ,t) + ^G~ 2 (k, k', t). (94) 

By using the identity GG^ 1 = 1, it follows that 

tfG-^kXt) = -G(k)~ 1 5G(k,k',t)G(k / ) _1 (95a) 

5G~ 2 (k,k',t) = 2G(k)~ 1 J 5G(k,k" ,t)G(k"y^G(k" ,k' ,t)G(k')- 1 . 

k" 

(95b) 

Upon promoting all quantities to the new variables (P,q), Hermitian con- 
jugation, SG(k, k',ty = SG*(k', k, t), obeys the mapping 

8G(k, k',t) -» <5G(q, P, t) ^ <5G*(k', k, i) -> £G*(q, -P, t), (96) 

which, by the definition of P given in (93a), follows from the interchange of 
k and k'. 

Undertaking a harmonic analysis, the mean grand canonical Hamiltonian, 
K = {H — fJ,N), must be expanded to second-order in the small deviations 
of (92), 

K = K<® + 6KSV + + ]^5Q} 03 SQ, (97) 

where the vector 5Q and its canonical conjugate 5V are given as 

ifit^P,*) = [SE(q,-P,t),57r(-P,t),Scj(-P,t),Su{-P,t)] (98a) 

<W + (q,P,t) = [5G(q,-P,t),5i;(-P,t),5<p(-P,t),5x(-P,t)}. (98b) 

In expression (97), the zero-order constant term, K^°\ is the static, uniform 
piece that can be dropped as it merely represents a constant energy shift, 
having no effect on any ensuing dynamics. Also, the first variation, 5K, is 
zero at the stationary points given by (64), thus leaving a quadratic Hamil- 
tonian in 5Q and SV. Consequently, 21 and 23 can be interpreted as the 
mass and spring matrices whose product has eigenvalues that determine the 
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low-lying excitations of the system. 

The second-order expansion of Ti — /mN" has been carried out explicitly in 
the appendix, which finds the 21 and 23 matrices to be 



a(q,q',P,t) 
<8(q,q',P,i) 



s M (q,P)£(q-q') c A/ (q,P) 



T 

Z M 



q'P) 



A(P) 



Sif(q,P)£(q - q') c^(q,P) 



(q'.P) 



BCP) 



(99a) 
(99b) 



where A(P) and -B(P) are 3x3 matrices [compare (A. 23)]. Finally, along 
with their corresponding transposes, and c^-, the off-diagonal terms are 

c M (q,P) = [2t G v(q,P) Oo^(q,P) ^MqP)] (100a) 
c K (q,P)=[S & (q,P) 9Wq,P) » s „(q,P)]. (100b) 



8.2 Quasi-Boson Interpretation 

As suggested by (99), it is natural to separate the diagonal or noninteracting 
part of Krpa from the off-diagonal piece so that 12 



Krpa = Kq + Ki n t, 



with 



K int = ~6V* 



SM 

A 



o a 



(101) 



(102a) 



For the noninteracting case involving only Kq, it is simpler to introduce the 



The subscript denotes the diagonal part of the matrix. 



34 QUANTUM INSTABILITY OF ATTRACTIVE BOSE SYSTEMS 
multiplicative canonical transformation 



/5Ti\ 




fSE\ 








5ir 




5tt 








5a> 




5lo 








\Su) 




V") 








[6G\ 




fSG\ 






\ 


5^ 




5ip 






5(j) 




5(j) 








W) 










J 



(103a) 



(103b) 



resulting in 
K = - (8Z* 5tt* 5uj* 5u*)' 



/0| 












fSE\ 















5tt 








w| 







5lo 


V° 








2 




\5u) 



+ l\5Vf 



(104) 



Associated with each of the respective pairs, (5Yl,5G), (5ir,5ij)), (Scj, 
(5v, 5x), is a frequency given by 



n 2 (q,P) = V^/(q,P)s^(q,P) 

W *(P) = ^(P 2 - M ) 2 -^ 2 

= e — 1\x 
tux = £ — 



(105a) 

(105b) 
(105c) 
(105d) 

An explicit calculation of 2 yields an interesting interpretation. From 
(A.25) and (A.35) we have 

SM (q, p) = |Gq: 1 G:V+ + (106a) 

s^(q,P) = 4G+GL(w+ +£«;_), (106b) 
where, in addition to the form of G in (64a), we have used the definition 13 

w± = w*(q±) = ^(4-/i) 2 -^ 2 . (107) 
Multiplication of sm and sk obtains 

^SMSK = LO++L0-. (108) 
Recalling the definitions of q and P in (93) leads to 

2 (q,P) = cu(k)+cu(k'), 



(109) 



3 The ± subscripts denote coordinate shifts of ±P/2 as in f± — /(q± P/2). 
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thus suggesting that u> and f2 2 are one an d two free quasi-boson energies, 
respectively. Although still incomplete, it is nonetheless possible to ex- 
trapolate these results to the more general case. Because they arise from 
an effective Lippmann-Schwinger equation, the eigenfrequencies of the full 
problem are identified with the energy of two interacting quasi-bosons. Ac- 
cordingly, discrete, real-valued frequencies represent bound states whereas 
a continuous range represents the scattering continuum. 



8.3 Eigenspectrum of 21 

Before trying to solve for the spectrum of 21- 25, it is helpful to first consider 
the solution for the eigenvalues of 21 alone, for this simplified case serves to 
illustrate the approach in the more general problem. Furthermore, a zero 
mode of either 21 or OS implies its presence in their product. For instance, 
consider a discrete basis where 21 is diagonal. If 21 has a zero eigenvalue, then 
multiplication by 25 reveals a product having an all-zero row, thus proving 
the persistence of the zero mode. Therefore, the analysis of 21 can also serve 
as a check on the results obtained thus far since the Goldstone mode is built 
into the equations from the particle number symmetry of the Hamiltonian. 

Throughout, our attention is primarily concerned with the eigenspectrum 
for the tp = solution. We begin with the eigenvalue equation 14 

21<5P = n a 5P, (110) 

which, with the aid of (98b), (99a) and (100a), can be written as four linearly 
independent equations in the variables 5G, dtp, 8(f), and 5\'- 

J s M (q)S(q " q')^(q') + 2l G ^(q)^ + 2l Gx (q)5 X = ^a<5G(q) (111a) 

q' 

i&w6i/> = n a 5i{) (111b) 

2l 0G (q')<5G(q , ) + % H 5(j) = n a S4> (111c) 
j 2l xG (q')<5G(q / ) + K XX 5 X = n a 5 x . (Hid) 

q' 

In addition to the notational conciseness of dropping the P indices, we have 
used (A. 30) and (A. 31) in setting 21^ and 2l^, x to zero for the tp = case. 



14 The subscript a is used to avoid confusion with the eigenvalue label used for 21 • 23 
(see Sec. 8.4). 
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Elimination of the fluctuations 5(f) and 5\ leads to a single equation in <5G(q), 
J {[ SM (q) - QJ 5(q - q') + C(O a )F(q, q')} ^(q') = 0. (112) 

q' 

Using the expressions for 21^, 2l xx , Sto^ and 21g x i n A. 6, V(q, q') simplifies 
to a single separable term, 

%q')=y(q)2/(q')> (H3) 

where we have defined 

(^'nMw^n^TT^ (114> 
»(q) = 5 (jfff'c: 1 + iY (115) 

Identifying the operator in (112) as 0(q, q') = [s^(q) — r2 a ]5(q — q') + 
£(fi a )y(q, q'), its spectrum is found by solving for the values of Q a at which 
an inverse, O -1 , does not exist. Without loss of generality, it may be as- 
sumed that 0~ l has the form 

o-V, q') = -, k o W " <*') + W> i'), (H6) 
sm (q ) - o 

where <S(q",q') is obtained from the requirement J C(q, q")C~ 1 (q // , q') = 
5(q — q'). Performing the necessary multiplications, this condition is met 
only if 

[s M (d)-n a ] 5(q, q ') + C(na) / ^(q,q")5(q",q') 

q " (117) 
= -C(O a )y(q,qO 1 

«M(q J - O a 

Defining a new matrix 

T(q, q') = [s M (q) - Q a ] S(q, q') [s M (q') ~ ««] , (118) 
gives the Lippmann-Schwinger equation used in scattering theory: 

r(q,q')+C(n«) / %q")^ rir T(q",q') = -C(^)^(q,q')- 

q" 

(119) 

Because this condition must be satisfied in order to invert the operator in 
(112), values of Q a are sought for which (119) has no solution. According 
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to the Predholm alternative [25], (119) cannot be solved if there exists a 
solution to its homogeneous counterpart, 



— L- T(q,q')+ / V(q,q") ± r(q",q')=0. 



(120) 



q" 

Thus, at least part of the spectrum is found by the {O a } that satisfy (120). 

As was done in the two-body case, we take a separable form for the T- 
matrix, 

T(q,q') =ty(q)y(q'), (121) 
which upon substitution into (120) gives 

1 + / " { f = 0. (122) 



C(n o ) 7 SAf(q) - 

q 

An exhaustive solution for the O a 's is not attempted here as the non-zero 
elements of the spectrum reveal little in acquiring the eigenfrequencies of 
21-03. However, a zero mode of 21 does prove its presence in the product, 
prompting us to check whether (122) may be solved for Q a = 0. Since the 
Goldstone mode occurs at zero total momentum, the required expressions 
are obtained from the P = form of (115) and (A. 25), thus giving 

y(q,P = 0) 2 _ 1 f {q 2 -vf 



s M (q,P = 0) 2 
q q 

_d_ 



3/2 



C(o) 



(123) 



Obtaining the second line requires the explicit form of J D(q) in (66c), 
whereas the last line is derived from a combination of (64c) and (114), thus 
verifying that the Goldstone mode is present in 21. 

In addition to Q a = 0, the inverse matrix, O" 1 , will not exist if 

= S A /(q,P) 



= \ rHr 7m? 7t( u + (124) 

Due to its dependence on q, this solution represents the continuum. More 
importantly, the system's excitations can be found from the eigenvalues of 
21 • 53 using this same approach. 
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8.4 The 21 • 53 Product 



For the full excitation spectrum, it is necessary to first obtain the matrix 
product which, by using the definitions (99), is 15 



2l-Q5(q,q / ,P)=M(q,q / ,P) 

"sA/(q)sAr(q')<Hq - q') + c M (q) c^(q') sjvf(q) c K (q) + c M (q)-B ' 

^r(q') *Jsr(q') + ^&(q') / &(q")c*: (q") + .4 • J3 



M GG (q,q',P) M G ^(q,P) M G0 (q,P) M Gx (q,P)' 

M^ G (q',P) M # (P) M # (P) M^ X (P) 

M 0G (q',P) M^(P) M 00 (P) M 0X (P) 

M xG (q',P) M^(P) M X0 (P) 



(125) 



Through (64a), (100), (106), (107) and (115), the matrix elements are ex- 
plicitly given by 

M GG (q, q') = (u+ + u;_) 2 <5(q - q') + (A 2 a? + g 2 )y{^)(G' + + G'_) (126a) 



M G ^(q) 

^Gx(q) 
M 0G (q') 



-(A 2 a 2 + <7 2 Mq)^ 

+ G= i)( w+ + u _) _ Aai y( q) ( e _ 2 ^) 

+ GI 1 )^ + w _) - 5 y(q)( £ - 2 M ) 



-(A 2 a 2 + 5 2 )(G / + + G'„)^ 
-4Aa 1 y(q , )G' + G , „K+^) 
-Aai (G" + + G'_)(e-2 M ) 
M xG (q') = -4 5 y(q')G' / + G / _(a;; W_) 
- 5 (G / + + G'„)(e-2^) 
= (P 2 - fi) 2 -e + (A 2 a 2 + g 2 )^ 2 
M n = Aai(P 2 -3/i + e + e)^ 
M^ x = 5 (P 2 - 3/i + e + e)V 



(126b) 
(126c) 
(126d) 
(126e) 

(126f) 

(126g) 

(126h) 
(126i) 
(126j) 



15 It is understood that all quantities within the matrix depend on the total momentum 
P = P'. 
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= Xai(P 2 - 3/i — £ + e)ip (126k) 

M H = X 2 aj (2h + V 2 ) + (e - 2fi) 2 (1261) 

M 0X = M x4> = Xaig (2ii + V> 2 ) (126m) 

M xi , = g(P 2 - 3/i - £ + e)il) (126n) 

M xx = <? 2 (27 1 + V 2 ) + (e - 2/x) 2 , (126o) 
where we have conveniently defined 

h^\j i/(q)(G+ + GL). (127) 
q 



These elements are used in the set of equations corresponding to the eigen- 
value problem, M 5P = Q 2 5P: 

[ M GG (q, q") 5G(q") + M G ^(q) + M G0 (q) 5<P + M Gx (q) 5 X 
i> (128a) 

= fi 2 5G(q) 

J M^ G (q") <5G(q") + M # 5V + 5<j> + M^ x <5x = ^ (128b) 

q" 

/" M^ G (q") 5G(q") + <fy + M H 54> + M 0X <5 X = fi 2 <ty (128c) 

q" 

y M xG {q") SG(q") + M x ^ <ty + M x0 ^ + M xx 5 X = ^ 2 S X - (128d) 

q" 

Up to this point the equations are general, describing the excitation spec- 
trum of any stationary solution. For the tp = case, Mq^, M^q, M^, 
M^, M^p x and M x ^, all vanish, reducing the above system to only three in- 
dependent equations. Further, we restrict the discussion to long wavelength 
excitations, P — ► 0, assuming the general case to have similar properties. 
With these simplifications, the molecular field displacements are obtained 
in terms of 5G: 

W = \J [(^ " M xx ) M^ G (q") + M^ x M xG (q")] <JG(q") (129a) 

q" 

6 X = j f - M H ) M xG (q") + M x4> M 0G (q")] <JG(q"), (129b) 
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where after defining 



ft 2 = n 2 - (e - 2fi) 2 

n 2 = n 2 - (e - 2fif, 

the denominator is expressed as 

d = (O 2 - M H ){<tf -M xx )- MfoM x<f , 

= n 2 n 2 - 2 g 2 h n 2 - 2\ 2 alh n 2 . 



(130a) 
(130b) 



(131) 



Substitution of (129) into (128a) obtains an effective matrix equation for 



SG, 



J [HcO J(q - q") + y(q) r(q") + *(q) S (q")] *G(q") = 0. (132) 



In addition to y(q) given in (115), the terms appearing in this expression 
are identified by the following: 



h(q) = 4u;(q) 2 -fi 2 

z(q) = G(c l )- 1 u(c l ) = 2(q 2 - P L-0 



(133a) 
(133b) 



r(q") = 2 \\ 2 a 2 + g 2 + ~ [\ 2 a 2 (e - 2 M ) 2 Q 



r>2 



+ 5 2 ( £ - 2„) 2 tt 2 - 2A 2 a 2 5 2 /i (6 - e ) a \ G(q" 



+ ^ [A 2 a 2 (6 - 2/i) fi 2 + g 2 (e - 2 M ) O 2 ] y(q")G(q")Vq") 



(133c) 



s(q") = [A 2 a 2 (6 - 2^) tt 2 + g 2 (e - 2/x) ft 2 ] G(q") 



+ ^(A 2 a 2 2 + < 7 2 2 ) y (q")G(q")^(q / '). 



(133d) 



The problem to be solved is to find all such that Eq. (132) has no solution 
except the trivial case 5G = 0. Along with a discrete set of eigenvalues, 
a continuum arises simply from the structure of the inverse matrix, as was 
indicated in the analogous calculation for 21. 
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8.5 Point Spectrum at Long Wavelength (P = 0) 

Proceeding in analogy with the method used for 21, the operator in (132) 
can be written as 

0(q,q") = MqH(q-q") + [y(q) 4^ 



s(q") 



(134) 



If it exists, the inverse, O 1 , has the form 



r(q') 

*(q') 



h(d>) ' 

(135) 

where T is some 2x2 matrix, yet to be determined. From the identity 
Jq,, 0(q, q") O-^q", q') = <5(q - q'), we have 

r(q"Mq") f r(q")s(q")- 



T+ 1 + 



/i(q") 7 /i(q") 
q" 

a(q")l/(q") /■ s(q")*(q") 



T = 0. 



(136) 



q" q" 

By the discussion following Eq. (119), there can be no solution if 

r(q)y(q) f s(q)z(q) 



Hi) 



l 



q q 

r(<i)z(q) f s(q)y(q) 



h(d) 



h(q) 



0. 



(137a) 



(137b) 



Hence, the values of O that solve (137) determine the point spectrum of 
21-03. 



8.5.1 Verification of the Goldstone Mode 

Because of the particle number symmetry in the original Hamiltonian, a 
zero-frequency mode must emerge as a general feature of this description, 
independent of specific parameter values. For a proof, it is easiest to take 
= then see whether Eqs. (137) are satisfied. Utilizing the variational 
condition (64c) along with the appropriate partial derivatives of (66c), all 
required integrals are expressed in terms of //, e, X 2 a\ and g 2 : 

G(q)y(q) 1 d 



h(q, O = 0)" 4£ dfi 
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= _l_ g 2 (e-2^ + X 2 a 2 (e-2^ 

2[ 9 2 (e _ 2 ^ )+A 2 a 2 ( e _ 2 ^)] 2 



G(q)*(q) 1 



DM 



/i(q, = 0) 2£ 
q 

_ 1 (e - 2/x) (e - 2//) 



2 5 2 (e - 2/x) + \ 2 a\ (e - 2/x) 

j/(q)G(q)MqMq) 1_0 /" 

fc(q, = 0) 8 9£ 7 

q 

1 (e - 2/i) (e - 2fi) 



8g 2 (e- 2fi) + A 2 a? (e - 2/x) 
y(q)G(q)Mq)z(q) 1 



/i(q,O = 0) 4 



(138b) 



(138c) 



Ji. (138d) 



Also, the = forms of d, r and s are 

d(fi = 0) = (e - 2^) 2 (e - 2,u) 2 + 2g 2 Ji (e - 2^i) 2 

+ 2A 2 a 2 /i (e - 2^f 
r (q, = 0) = ^ [ 5 2 (e - 2 M ) + A 2 a 2 (e - 2/x)] 2 G(q) 

- ^ (e - 2/i) (e - 2/i) [g 2 (e - 2/x) 

+ A 2 a 2 (e-2/x)] y(q)G(q) 2 u,(q) 



(139) 



(140) 



s(q, = 0) = -- d (e - 2fj,) (e - 2 M ) [ 5 2 (e - 2 M ) + A 2 a 2 (e - 2//)] G(q) 



g 2 (e - 2 M ) 2 + X 2 a 2 (e - 2^) 2 j y(q)G(q) 2 W (q). 

(141) 



Finally, a combination of (138)-(141) satisfies the eigenvalue conditions (137), 
confirming the presence of the Goldstone mode for arbitrary parameter val- 
ues. 
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8.5.2 Eigenfrequencies in the Zero-Range Limit 



We now solve for the nonzero solutions to the eigenvalue equations in the 
zero-range limit, A — > 0~. First, r and s are expanded to the appropriate 
order in A, then used in (137). Defining I\ = k\~ 3 + const., all subsequent 
expressions may be conveniently written in terms of k. Accordingly, the 
denominator has the limit 

1 £ 2 A + 2^A 2 + ^^^A 2 T 



A 5 



d a^o- 



1 



2kq 2 



An 2 af 



When substituted into (133c) and (133d) this leads to 

K — 1 



r(q) 



A^O" 



2 9 
2 £ 



a 



4 (e - 2/i) A 



KOI 



G(q) 



(142) 



+ 



4 r0 2 A 2_l (e _ 2/i) A 3 



(143a) 



+ 



2^,4 



y(q)G(q) 2 o;(q) 



s(q) 



A^O" 



1 4\3 
— fi 2 A 2 G(q) y(q)G(q) 2 u;(q). 



(143b) 



Next, it is necessary to calculate the small A versions of (138). However, in 
this case the variational condition cannot be utilized since Vt ^ in general. 
It is easiest to write these expansions in terms of the following integrals 16 



J 



1 = 7T 



q [q 



i 



n' 2 



1 



7T 

"4 ^y 2 



n> + y/-n + n> 



(144a) 



J: 



2 — 7: 



(q 2 -n?-e 



( q 2 -fiy-n' 2 



dq 



nyJW 2 



y/-H-M -y/-fl + n' 



(144b) 



16 These expressions are obtained by first factorizing the integrand denominator as 
(q + — Q')(q — i\/~-JT^ViJ){q + i^—j-i + Q')(q — i\/~-JTTTV), then choosing either 

the upper or lower semicircular contour. 



44 QUANTUM INSTABILITY OF ATTRACTIVE BOSE SYSTEMS 



where f2' 2 = £ 2 + Q 2 /4. To the required order in A, the analogs of (138) 
may now be expressed as 



f g(g)y 

J Hq 



G(q)iKq) ^ J_ 

) A^O- 167T 2 



Ji (145a) 



— TT^ — oT Tr h F~2 J 2 (145b) 

q 

y(q)C(q) 2 w(q)y(q) 1 



/i(q) A->0- 8A 

q 

y(q)C(q) 2 w(q)2;(q) K 



+ const. (145c) 

» — ^ + const.. (145d) 
h(q) x-^o- 4A 3 

q 

With Equations (143) and (145), the eigenvalue conditions (137) become 



r (qMq) J\ 2 , / 1/lfi n 

n^T^-s^ 20 ^- 1 (146a) 

r(q) ; (q) >e -2,-(^ + A)l o 2 = (146b) 

h(q) A^o- P V8vra b!7 4tt 2 / a 2 e V ; 



g(q)z/(q) 

/i(q) A-+o- 
s(q)z(q) 



/i(q) A^O- 



C(A 2 ) = (146c) 
-1 + 0(A), (146d) 



thus showing that the last two are automatically satisfied for zero range. 
Aside from the Goldstone mode, the other eigenfrequencies are obtained 
most easily by solving (146a) and (146b) for J\ and J2 then equating their 
product to that given by a direct calculation using Eqs. (144): 

7T 2 A A , si 47r 3 a? 1 

JlJ2 = __ = 32,* f( e -2 M ) f? -^-l ?f . (147) 
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From the definition of (130a) and the solution of (147), two distinct 
eigenfrequencies finally emerge as 

167ra 2 



(e - 2fi) z + 



a bg 



1± \ 1 



(e - 2/x) 



(148) 



depending on the density through fi only. In addition to a stable response 
frequency of the many-body collective, real values of Q may be interpreted 
as the binding energies of two interacting quasi-bosons. On the other hand, 
scattering states are implied by complex Q, with the imaginary part quanti- 
fying damping in the corresponding collective excitation. It turns out that 
for the case depicted in Fig. 6, both roots of Q are complex, having real 
parts that lie at ~ 23 neV, well within the scattering continuum. 



8.6 Continuous Spectrum of the Lower State 

Unlike the discrete elements, the continuum is obtained directly from the 
construction of the inverse. Specifically, the condition /i(q) = also pre- 
cludes 0~ l from existing, thus giving the continuous part of the spectrum, 

n = 2w(q) 

, (149) 

= 2vV-^) 2 -e 

Since the form of h(q) is the same regardless of whether ip = or ip ^ 0, 
Eq. (149) remains valid for the continuum of the entire lower state. At long 
wavelength (P = 0), these excitations are parametrized only by q 2 , with the 
lowest energy obtained for q = 0. Consequently, half of the lowest quasi- 
boson energy, 0/2| g= o, plus the energy per particle gives the continuum 
boundary: 

E exc = e + vV-£ 2 - (150) 
As shown in Fig. 6, this threshold lies above the tfj = solution (/x ^ £), 
but intersects with the collapsing state at the critical point. Thereafter, 
the boundary coincides with the ^ / piece (/i = £), indicating that the 
continuum contains the zero mode for the second piece of the collapsing 
solution. 

It is now possible to give the decay found in Section 6 a more complete 
physical explanation. Note that the real part of the energy per particle lies 
inside the continuum of the collapsing state. Due to energy conservation, 
this upper solution can only evolve into a state at the same level, suggesting 
the coherent decay be interpreted as a transition into collective phonon 
excitations which inherit the same collapsing behavior as their associated 
ground state. Hence, the condensate is initially lost through decay into 
the phonons of the lower level, as opposed to a physical loss of particles. 



46 QUANTUM INSTABILITY OF ATTRACTIVE BOSE SYSTEMS 

Nevertheless, as the system collapses, qusiparticles will emerge causing a 
kinematical atom loss akin to the semiclassical recombination models [5, 6]. 
However, the previous analysis is only valid for low-temperature systems, 
and therefore may not be applicable to these kinematic effects. 



Figure 6. A figure showing the same curves as in Fig. 5, 
where the hatched regions indicate the continuum of modes 
belonging to the collapsing lower state. As discussed in the 
text, the excited "false vacuum" always lies in the continuum 
of the lower state, thus indicating that the decay represents a 
transition into collective phonon excitations of the collapsing 
solution. Note that the density of states for tp ^ is greater 
than that for tp = 0, as shown by the hatched lines above each 
region. Additionally, the real parts of the discrete spectrum 
given by (150) lie relatively far up in the continuum at ~ 23 
neV. 



9 Summary 

We have considered the case of an atom-molecule condensate in which the 
interactions were attractive yet the effective scattering length was positive. 
This situation raised the question of whether the condensate collapsed due 
to the mutual attractions or remained stable in accordance with the positive 
scattering length. Starting with a two-body analysis, a separable potential 
was used to realistically model the interparticle interaction. Due to its suc- 
cess, this separable form was implemented in the many-body Hamiltonian. 
Equations of state were then obtained from the application of a variational 
principle that utilized a Gaussian trial wave functional for the many-body 
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state. Despite the positive scattering length, a collapsing solution was ob- 
tained, consisting of a relatively low-density piece having only a molecular 
component to its condensate. At higher density, there occurred a quan- 
tum phase transition after which the solution comprised both atomic and 
molecular condensate components. Only by allowing the chemical poten- 
tial to assume complex values could the experimentally observed case be 
obtained. As the phase of the order parameter, the chemical potential has 
an imaginary part quantifying a decay rate, assigned a physical meaning 
through a small oscillation analysis about the equilibrium solutions. Ex- 
panding around the stationary points revealed two discrete eigenfrequencies 
associated with the low-density molecular condensate solution. Moreover, 
the experimentally observed energy per particle lay within an excitation 
continuum of the collapsing two-piece lower state. From energy conserva- 
tion, the decay was interpreted as a coherent process corresponding to the 
evolution of the observed case into the excitations of the lower solution. 
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Appendix A Derivation of the RPA Hamiltonian 



In the appendix, we expand the expectation value of the Hamiltonian (57) 
to second-order in the variations given by (92), thus deriving the 21 and 23 
matrix elements used throughout Section 8. 



A.l Expansion of R, D, <3? , X 2 and S 



As appears in (57), the momentum space forms of the quantities are given 
by 



i2(k,k',t) 



1 



G _1 (k, k', t) + G(k, k', t) - S(k - k 1 ) 



+ 2 J S(k,k ,, ,t)G(k",k w ,t)S(k /,/ ,k / ,t) 



k",k" 



(A.la) 



D(k,k',t) = - 



-G-\k,k',t)-G(k,k',t) 



+ 2 J S(k,k ,/ ,t)G(k",k /,/ ,t)S(k /,/ ,k / ,t) 



k",k'' 



(A.lb) 



[E(k, k", t)G(k", k',t) + G(k, k", t)£(k", k', t)} 



1 

72 



^(k, t) = —= V(k, t) + zvr(k, t) 



(A.lc) 



The molecular field expressions are straightforward, following from (A.lc). 
Expanding K to second-order first requires an expansion of R, D, X 
and H in terms of the variations (92). Using (A.l) obtains 



R(k, k',t) = R.W (k, k',t) + R {1) (k, k', t) + R {2) (k, k', t), (A.2) 



where 



i?(°)(k,k',t) = i 



^(k^+GCk)-! 



<5(k - k') 



(A.3a) 



i?(k)5(k - k') 
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i?M(k,k',i) = - 



-G(k)- L G(k')- L + 1 



«5G(k,k / ,t) 
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(A.3b) 



i?( 2 )(k,k',t) = ^G(k)' 1 y 5G*(k,k ,/ ,t)G(k ,/ )~ 1 5G(k /, ,k',t)G(k / )~ 1 

k" 

+ 2^ 5S*(k,k",t)G(k")5S(k // ,k / ,t). (A.3c) 



Similarly, the first few orders of the other fluctuations are 



D(°)(k,k',t) = i 



^G(k)- 1 - G(k) 



<5(k - k') 



(A.4a) 



D(k)tf(k - k') 



5 (1) (k,k',i) 



fG(k)- 1 G(k')- 1 + 1 



<5G(k, k',t) 



(A.4b) 



- t [G(k) + G(k')] <5S(k,k',t) 
D^Oc^.t) = ^G(k)" 1 J 5G*(k 1 k",t)G{k")~ 1 5G(k",k',t)G(k'y 1 



k" 



+ 2^ 5S*(k,k // ,t)G(k")5S(k // ,k',t) (A.4c) 
k" 

i y [5Z* (k,k" ,t)5G{k" ,k' ,t) + 5G* {k,k" ,t)SZ(k" ,k' ,t)]. 



Expanding the products of mean fields gives 
^*^(°)(k,k',t) = ^ 2 5(k)5(k') 
tf**W(k,k , ,t) = -^<5(k') [<fy*(k,t) -iSir*(k,t)] 
+ I?p5(k) [5i>(k',t) + i5ir(k',t)] 



(A.5a) 
(A.5b) 



^*^ (2) (k,k',t) = - \5i/)*(k,t)6M,t) +i5ij*(k,t)5iT(k',t) 

2 (A.5c) 

- i#(k', t)5vr*(k, t) + for*(k, t)<j7r(k', t)] 



50 QUANTUM INSTABILITY OF ATTRACTIVE BOSE SYSTEMS 



**(°)(k, k',t) = \^\? <5(k)<5(k') (A.6a) 
(k, k',t) = -if) 6(k') [<fy(k, t) + i5n(k, t)] 



+ ^5(k) [5i/;(k',t)+i6<ir(k',t)] 

^ 2 \k,k',t) = - [5i>(k,t)dip(k f ,t) +i6ip(k,t)dTr(k',t) 
+ i8ip(k', t)5n(k, t) - 5ir(k, t)6ir(k', t)] 



+ ^<f>6(k) [^(k / ,t)+i5u;(k , ,t)] 

$*$( 2 )(k,k',£) = - [^*(k,t)^(k',t) + ^*(k,t)<5o;(k / ,t) 
- i84>{k!, t)8u*(k, t) + 5w*(k, t)6u(k',t) 



+ \xm [6x(k',t)+i5v(k',t)} 



X*xM(k,k',t) = \ [S X *(k,t)5 X (k',t) + i6 X *(k,t)6v(k',t) 



2 

- iS X (k', t)5v* (k, t) + <5z/"(k, t)8v(k', t)] 



(A.6b) 



(A.6c) 



$*$(°>(k,k',t) = ^cf> 2 5(k)5(k') (A.7a) 
$*$W(k,k',t) = ^4>6{k') [5<t>*(k,t) -i5uj*(k,t)} 



(A.7b) 



(A.7c) 



A*A(°)(k,k',t) = i x 2 ,5(k) ( 5(k') (A.8a) 
X*xW(k,k>,t) = ± X S(k>) [6 X *(k,t)-i6v*(k,t)] 



(A.8b) 



(A.8c) 



~^(k",t) = j=8(k») 

(A.9a) 

= -L(Aa 1 + 5 x)«5(k / ') 
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1 



= — [\ ai 5cl)(k",t)+g6xQt",t)] 

+ -5= [Aai5w(k", t) + ff^(k", t)] . 
v2 



(A.9b) 



Along with (A. la) and (A. lb), these expansions are used to calculate the 
second-order contribution from each term in (57). 



A. 2 Kinetic Contribution 

From Eqs. (A. 3), the second-order part of f(k 2 — /i)-R(k, k, t) is 

k 

J(k 2 -^)R^(k 1 k,t)= J [5G*(q,P,t)(A ; 2 -M) J R G2 (q ) P)JG(q J P,t) 

k q,P 

+ 5£*(q, P, t)(k 2 - //)i? E 2(q, P)«5£(q, P, t)] , 

(A.10) 

where the total and relative momenta are identified in (93). Furthermore, 
we have used the transformation 6G*(k', k, t) — > 5G*(q, P, t), with the sym- 
metry under inversion P — > — P incorporated by the following definitions 



(fc 2 -/^ G2 (q,q',P) 



(fc 2 -/i)i2 E2 (q,q',P) 



1 

16 



(q 2 , - ^GfGZ 1 
+ (q 2 -^GZ^J^q-q') 
(q 2 - M )G_ + (q 2 - M )G + lff(q-q') 



(A.lla) 
(A.llb) 



For notational simplicity, the ± subscripts denote a shift of ±P/2 in the 
argument of any function / such that 



/ ± = /(i±ip), /i = /(o'4p')- 



(A.12) 
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By (A. 5c), the second-order contribution from f (k 2 — fi)^!*(k, i)^(k, t) is 

k 

J (k 2 - (k,\c,t) =j(k 2 -^ [<ty>(k,i) 2 +67t(k,t) 2 ] 

[5r(p,t)(k 2 -/x)***^(p)^(p,t) 

+ <5vr*(P, t)(& 2 - p)***** (P)5vr(P, t)] . 

(A.13) 

In addition to using the fact that the fluctuations are real, we have defined 



(k 2 - V(P) = -(P 2 - M ) (A.14a) 
(A; 2 -^)^^ 2 (P) = i(P 2 - / u). (A.14b) 



A. 3 Molecular Contribution 

The second-order contributions from the molecular detunings are given by 
(e-2fi) J <S>*& 2 \k,k 1 t) + (e-2 f i) J X*X® (k, k, t) 

k k 

= j [^*(p,t)$*^ 2 (P)^(p,t) + < yw*(p,t)**$ w3 (P)(yw(p,t) 
p 

+ < 5 x *(P,t)A*A x2 (P),5 x (P,t) + fo'(P,On„ 2 (P)fo(P,t)]. (A.15) 
Reading off the results from (A. 7c) and (A. 8c), we identify 

$*^.(P) = fc*$ wa (P) = ie-^ (A.16a) 

X*X X 2 (P) = X*X V 2 (P) = -e - n. (A.16b) 
A. 4 Contribution from the Coupling 

The last integral of (57) gives the energy due to the coupling of atoms to 
both of the molecular states. In this contribution, there is an H*D term 
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T2 J a < k "- k + k '"v 

k,k',k" 



k + k' 



<°)(k",t) J D( 2 )(k,k , ,t) 



+ E* W (k",t)D ( - 1) (k,k',t) +H. c. 



/(k)£ 



k,k" 



<5G*(k, k"\ i) -G(k)- 2 G(k" / )" 1 <5G(k w , k, t) 
8 

+ <5£*(k, k'", t)2G(k"')5E(k"' , k, t) 
k + k' 



+ / <Hk"-k + k')/(^) 



k,k',k" 



^i^*(k",t) 



<5x*(k",i) 



i G(k)- 1 G(k')" 1 + 1 



G(k)- 1 G(k')" 1 + 1 



5G(k, k',t) 



5G{k,k',t) 



- Xai5uj*{k",t) [G(k) + G(k')] *E(k, k', t) 
-gSu*(k",t) [G(k) + G(k')] £E(k,k',t) j, (A.17) 



which can be more compactly written as 



q,q',P 



+ 



1 [<5G* (q, P, t)D G2 (q, q', P)5G(q, P, t) 

+ 5£*(q, P, t)D^ (q, q', P)<JS(q, P, t) 
'(P, t)I>G*(q, P)5G(q, P, t) + 5 X *(P, 0^ x (q, P, t)6G(q, P, *) 



q,P 



+ <^*(P, ^^(q, P)<5E(q, P, t) + 5v*(P, t)Dxu(<i, P)^(q, P, t) 



(A.18) 
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where 



£> G 2(q,q',P) =£ 



^(q,P) 



2/16 



/[q+^Pl-G^GZ 1 



1_\ 1 



-2 ^(-1 



(A.19a) 



<5(q-q') 



L> s2 (q,q',P) = £ 



+/[q_-P)G + 



^/(q)(^ G - + l 



%-q') 



£ Gx (q,P) = -|/(q) Qg^G^ + i) 

^(q,P) = -Aai/(q)(G+ + GL) 
D a ,(q,P) = -<7/(q)(G + + GL). 



(A.19b) 

(A.19c) 

(A.19d) 

(A.19e) 
(A.19f) 



Likewise, the second-order part from the S*^^ term is 



k,k',k" 



<°\k",t)^ 2 \k, k',t) 



+ 3*M(k",*)** (1) ( k , k ' ) *) +H. c. 

= | /(k)|[#(M 2 -^(M 2 ] 

k 

+ i J 5(k"-k + k')/(^^){[A ai ^(k",t)+^ x (k",f)] 

k,k',k" 

x [#(k')#(k,t) +#(k)^(k',t)] 

+ [\ ai 5u)(k",t) +gSu(k",t)) [#(k')(57r(k,t) +#(k)57r(k',t)]}. 

(A.20) 
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After multiplying out the integrand, this becomes 

J [^*(P,t)¥ 2 ^(P)^(P,t) +5n*(P,t)*l 2 (P)5TT(P,t) 
p 

+ 50* (P, P)^(P, t) + 5 X *(P, t)^ 2 xV ,(P)^(P, t) (A . 21) 

W(P,t)1'L(P)^(P,t) +<^*(P,t)*L(P)&r(P,t)" , 

where 

^ 2 (P) = -^ 2 (P) = i/(P) (A.22a) 
* 2 ^(P) = *L(P) = Aqi^/ Qp) (A.22b) 
* 2 xV ,(P) = * 2 W (P) = #/ Qp) . (A.22c) 



A. 5 Matrix Form of the RPA 



Expressed in matrix form, the grand canonical RPA Hamiltonian corre- 
sponding to Eq. (57) can be written 17 



Krpa 



^Q f (q) 



S^G(q') 
axG(q') 



^s(q') 
® W E(q') 
»^(q') 



2Mq) 


2l G 0(q) 


2tGx(q)" 










21^ 


2W 




%c<t> 


^xx . 


®E^(q) 


93s* (q) 


»E„(q)" 






237rf 


23 




23 aw 


®I/7T 







with the coordinates, 5Q, and momenta, defined as 
(JQt(q )P ) = [«5S(q,P) <frr(P) <Jw(P) <5i/(P)] 

5^(q,P) = [<5G(q,P) <ty(P) ty(P) *x(P)]. 



*e(q') 5 

(A.23) 

(A.24a) 
(A.24b) 



7 The dependence on the total momentum P has been omitted for notational clarity. 
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To find explicit expressions for the matrix elements, we combine the results 
of (A.ll), (A.14), (A.16), (A.19) and (A.22). 



A. 6 21 Matrix Elements 

As in the discussion following (64), we take /(k) = 1 with a cutoff at |k| = 
An 2 /b. Subsequently, the matrix elements are simply matched with the 
expansions, leading to the following identifications: 

2l GG (q, q', P) = 2R G 2 (q, q', P) + 2D G 2 (q, q', P) 

= I [(<£ - fi + OGfGZ 1 + (<£ - m + OGZ 2 G^] 5(q - q') 
= s A .f(q,P)5(q- q). 

(A.25) 

In the above case, for instance, the terms were obtained from (A. 11a) and 
(A. 19a). Similarly, the remaining elements are: 

2l G ^(q, P) = 2V; (q, P) = (A.26) 



2l G0 (q, P) = 2l^ G (q, P) = D g ^cl, P) 

21 Gx (q, P) = 2t xG (q, P) = D Gx (q, P) 



(A.27) 



2l^(P) = 2(k 2 - ^)***^(P) + 2* 2 , 2 (P) 



2l^(P)=2t^(P) = ^(P) 

2t^(P)=2l xV ,(P) = ^(P) 

= # 

2t^(P) = 2«D*V( P ) 
= e-2fi 



(A.29) 
(A.30) 
(A.31) 
(A.32) 
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a^(P) = 2y (P) = (A.33) 

2L y (P) = 2X*X y2 (P) 

/ (A.34) 

A.7 *B Matrix Elements 

®EE(q, q', P) = 2i? E2 (q, q', P) + 2D & (q, q', P) 

= 2 [(<£ _ M + £) GL + (<£ - ,x + £) G+] 5(q - q') (A.35) 
= ^(q,P)5(q-q'). 

®&r(q, P) = ®^s(q, P) = (A.36) 

^(q, P) = 23 w s(q, P) = £>£u,(q, P) 

= -Aai (G+ + G_) 



®E„(q,P) = ®„s(q,P) = £>s,(q,P) 

= - 5 (G+ + G_ 



»™(P) = 2(fc 2 - M )***^(P) + 2^ 2 (P) 



P 2 ~ M - £ 



2WP) = 2WP) = *L(P) 

= \a\ij} 

<B^(P) = » W (P) = tf 2 w (P) 

= # 

«8 ww (P) = 2*** wa (P) 
= e-2n 



(A.37) 
(A.38) 
(A.39) 
(A.40) 
(A.41) 
(A.42) 



» ul ,(P) =» W (P) =0 (A.43) 
® W (P)=2X*^(P) 

(A.44) 

= e - 2/i. 
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